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Abstract 

This simulation study investigates the dependence of the structure of dry Nafion®-hke ionomers 
on the electrostatic interactions between the components of the molecules. In order to speed 
equilibration, a procedure was adopted which involved detaching the side chains from the backbone 
and cutting the backbone into segments, and then reassembling the macromolecule by means of 
a strong imposed attractive force between the cut ends of the backbone, and between the non- 
ionic ends of the side chains and the midpoints of the backbone segments. Parameters varied 
in this study include the dielectric constant, the free volume, side-chain length, and strength of 
head- group interactions. A series of coarse-grained mesoscale simulations shows the morphlogy to 
depend sensitively on the ratio of the strength of the dipole-dipole interactions between the side- 
chain acidic end groups to the strength of the other electrostatic components of the Hamiltonian. 
Examples of the two differing morphologies proposed by Gierke and by Gebel emerge from our 
simulations. 
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I. INTRODUCTION 



Typical proton-conducting polymer electrolyte membranes (PEM) for fuel cell applica- 
tions consist of a perfiuorinated polymeric backbone and side chains that are terminated by 
strongly acidic sulfonic head groups. They provide the transport medium for the protons 
generated during the anodic oxidation reaction of the fuel. PEM fuel cells offer numerous 
advantages as power generators for emission-free vehicles, portable applications, and indi- 
vidual homesi. When such ionomers are exposed to water or humid air, the acidic groups 
terminating the side chains dissociate into SO3 groups on the chains, and protons in the 
aqueous sub-phase. The sub-phase extends through the membrane, creating a bi-continuous 
nano-phase separated network of aqueous pores and polymer. The percolation structure 
and its ability to support proton diffusion have to meet strict requirements for a PEM to 
be considered as an 'ideal' membrane for fuel cells: a) It must form a chemically, mechani- 
cally, thermally and morphologically stable structure over a wide range of water content, b) 
show steady performance and tolerance of elevated temperatures up to 150°C at high proton 
conductance, and c) be impermeable to gases and contaminating ions. Development of new 
PEM membranes that could meet all these demands is a primary task for many research 
groups involved in fuel-cell studies. 

One extensively used PEM material is Nafion®, a perfiuorinated polymer from DuPont. 
The integrity and structural stability of Nafion membranes is provided by the polytetra- 
fiuoroethylene (PTFE) backbones, but it is the hydrophilic clusters of side chain material 
that facilitate the transport of ions in the membrane. The microstructure of these clusters 
consists of an inter-facial region of solvated perfiuoroether side chains separating the polymer 
matrix from the more bulk-like water found in the pores. In the presence of water, one finds 
a nano-separation of hydrophobic and hydrophilic domains. According to the early work of 
Gierke et al?'^ these pores are inverted micellar spheres with diameters around 4 nm, and 
are interconnected by channels with a diameter and length of about 1 nm. 

Recent experimental investigations suggest that the spherical shape and uniform spacing 
of these pores are serious oversimplifications^. New small-angle X-ray scattering (SAXS) and 
small-angle neutron scattering (SANS) experiments by Gebel et al.^ for Nafion membranes 
in different swelling states have given a different picture for the morphology. Their results 
suggest that elongated bundles of hydrophobic backbone material with diameters of about 4 
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nm and lengths greater than 100 nm are surrounded by sulfonate groups and water molecules. 
In other words, hydrophobic material is surrounded by hydrophilic domains. 

The fact that controversy continues to surround the nature of the hydrated morphology 
and the shape of the domains suggests that there is still room for theoretical investigation 
and computer simulations to address this issue. Ab initio electronic structure calculations 
of polymeric fragments in contact with water- and quantum molecular dynamics studies on 
model PEM systems^ have provided a basis for understanding the energetically favorable 
conformations of Nafion chains and the molecular ingredients of the conduction process^*^. 
For systems large enough to contain several examples of hydrophilic inclusions and their 
interconnecting networks, however, treatment of the polymer in an ab initio manner, i.e. a 
full electronic treatment with molecular orbital theory, is not computationally feasible. 

The objective of this paper is to develop a deeper understanding of the driving mechanism 
behind the nanophase aggregation of sulfonate head groups, and to test the extent to which 
coarse-graining methods can be applied. We attempt this by using a simple dipolar head 
group model for ionomer side chains. We show that the dipole-dipole correlations strongly 
depend on the molar fraction of acidic groups and dielectric permittivity of the membrane. 
Globular micelles formed in the case of zero dipole moment evolve to percolated cluster 
structures for a non-zero dipole moment of head groups when an artificial attraction between 
the end groups of side chains is introduced. Accounting for the dipole moment of sulfonate 
head groups thus changes the conformational structure of clusters from a Gierke model of 
spherical aggregations to a Gebel model of elongated aggregations. The determination of the 
parameters necessary for success in modeling nanophase aggregation will help us understand 
the essential role played by the electrostatic interactions of the head group in determining 
microstructure formation in PEM materials. 

The rest of this paper is organized as follows. After a brief review of segregation models 
in section II, we present our dipolar model for side chains in section III. Section IV is devoted 
to our simulation results for dry PEM ionomers. We conclude with a discussion of results 
and planned future work in section V. 
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II. NANOPHASE SEGREGATION MODELS 



Significant experimental and theoretical efforts have been put into the characterization of 
the microstructure of dry and swollen Nation membranes. Although the fact of microphase 
segregation and its influence on the thermodynamic and transport properties of Nation mem- 
branes are well established, some uncertainty regarding the morphology of the subphases 
in swollen Nation remains. It has been generally assumed that the hydrophilic phase is 
continuous, in accord with the experimental evidence of high conductivity and water perme- 
ability. The conventional model of microphase segregation in Nation membranes, put forward 
by Gierke et al^'^, where the protons and water diffuse between globular clusters through 
cylindrical narrow channels, was further modified by Mauritz and Rogers^^ and Eikerling et 
alr^. Also, loselevich et ai— developed a lattice based micelle-channel model to describe a 
connecting structure of a water bridges between the cages of the micelles. However, there is 
as yet no direct experimental evidence for the existence of the channels connecting the clus- 
ters in Gierke's model. Other models explaining the properties of hydrated Nation exploit 
bilayer, lamellar, and sandwich-like structures for nanophase separation in PEM membranes. 
Gebel^ describes the Nation membrane as an aggregation of polymeric chains forming elon- 
gated objects (simplified as cylinders), embedded in a continuous ionic medium. At larger 
scales, those aggregates form bundles with definite orientational order. This new picture of 
multi-scale structure can explain membrane swelling as a continuous process from the dry 
state to a colloidal suspension. 

Scanning calorimetry experiments on cast Nation filmsi^ led to the conclusion that a 
continuous water phase exists only when the relative humidity of the membrane exceeds 90%. 
This conclusion contradicts the simulation results of Khalatur et al.—, where formation of 
continuous channels at very low hydration, and even in almost dry membranes, is reported. 
This is certainly consistent with the conductivity measurements in Nation films with very 
little water present^. In contrast to all the above-mentioned network models, Vishnyakov et 
ai— showed that the water clusters do not form a continuous hydrophilic subphase. Instead, 
the cluster-size distribution evolves in time as a consequence of the formation and rupture 
of temporary bridges between clusters. 

Molecular dynamics (MD) studies, using both classical force fields and an ab initio ap- 
proach, appear to support the idea of irregularly shaped ionic clusters of sulfonic groups 
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connected by smaller channels of similar composition embedded in the hydrophobic cova- 
lent matrix. However, the characteristic length scales of the morphological conformations 
suggested by experimental data are currently inaccessible to atomistic modeling. To probe 
the morphology at the nanometer scale, a mesoscale model is required. 

III. DIPOLE MODEL FOR SIDE CHAINS. 

In this work we are concentrating on the effects of electrostatic interactions among the 
elements of the side chains, and so in order to avoid the excessively slow processes resulting 
from chain entanglement we consider a segmented system in which the side chains are 
detached from the backbone. The later introduction of very strong interactions between the 
ends of the side chains and the backbone, which have the effect of knitting together the 
segments into a single macromolecule, do not seem to affect our results significantly. 

A generic acidic side chain of Nation is commonly denoted as 
CF30CF2CF(CF3)0(CF2)2S03H. In order to expedite these simulations, we adopt a 
simplified model for the side chain by using a united-atom presentation of the CF2 and 
CF3 groups^'i^'i^ and the SO3" sulfonate head groups^^. These groups are modeled as 
Lennard- Jones (LJ) monomers. We employ the same united-atom approach for backbone 
chain segments. 

In the framework of our dipole model all the partial charges on the ether oxygen, carbon 
and fluorine atoms of the side chain are set to zero. The fiuorocarbon groups of the backbone 
skeleton are also assigned zero partial charges. The electrostatic charge is located entirely 
on the sulfonate head group, SO3 , and on the hydrogen ion H"*" permanently bound to it. 
The dipole moment produced when the proton and the sulfonate group, modeled as LJ 
monomers, carry the full formal charges of +e and — e respectively, was reduced by a factor 
Dr in some of our simulations. Though our model is not applicable for the measurement of 
proton diffusion in PEM membranes, it is a good starting point for a step-by-step exploration 
of the nanophase morphology in dry PEM materials. 

The overall inter-monomer interaction in our model system is a combination of elec- 
trostatic Coulomb interactions between each pair of charged particles, 12-6 Lennard- Jones 
(LJ) potentials between neighboring monomers, and intrachain stretch and bend energies. 
The Lennard- Jones interaction parameters for the side chain are chosen to agree in most 
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instances with the Nafion model of Paddison22 . We make a further simplification for the 
LJ parameters of the side chain by assigning a universal value of a, chosen to be 3.5 A for 
all monomers. The total potential energy of the system is then given as 

^total = ^Lj + Eq + ^bond + ^angle' (1) 

where -Elj, Eq, E-^^^^, and E^^^^ are the Lennard- Jones, electrostatic, bond-stretching 
(two-body term), and angle-bending (three-body term) components of the total energy, and 
are given by 

ElAt) = AeL.M''/r'' - a'/r') (2) 
Eq = Y.^%^ (3) 

E^^^^{R) = \h{R-B^f (4) 

^angleW = ^^^(^-^o)' (5) 

We omit consideration of the dihedral potential term (four-body term), as this would be 
expected to make only a small contribution to the energetics of the single chain confor- 
mations in these dense polymeric systems. (This assumption was tested in sample sim- 
ulations in which this term was restored, and found to have little effect on our results.) 
The following parameters for the polymer chains were used: an equilibrium bending an- 
gle Oq = 110°, an equilibrium bond length Rq = 1.54A (0.44cr), a bending force constant 
kg = 120 kcal/mol deg^, and a stretching force constant kh = 700 kcal/mol A^. The last two 
constants correspond to the thermal energy fc^T at room temperature for a bond- length 
displacement of 0.1 A and bond-angle displacement of 1° from the equilibrium position. 
The bond length between the SO3 head-group monomer and the hydrogen ion was chosen 
to be cr, and hence 0.35 nm. The LJ interaction coefficient Elj was chosen to be SksT for 
hydrophobic-hydrophobic (HH) interactions, and IksT for both hydrophobic-hydrophilic 
(HP) and hydrophilic-hydrophilic (PP) interactions. In the latter case The LJ potential 
was modified to be purely repulsive by truncating at its minimum, where r = 1.1224cr, and 
raising it by the addition of an amount elj- The reason for this choice is that the interaction 
between protons is not well represented by van der Waals forces, which arise from virtual 
dipole-dipole interactions, and also that such attractive forces are small in comparison to 
the Coulomb interaction. The normalized charge Qi (/ = i,j) in Eq.(3) is +1 for protons 
and -1 for sulfonate groups. 



We used this coarse-grained model to perform molecular dynamics simulations of ensem- 
bles of ionomer segments. Because chain entanglement greatly slows the rate of equilibration 
in polymers, we adopted the strategy of temporarily cutting the system into more mobile 
segments, which were later reassembled. In order to increase the diffusional rate of our model 
ionomer we detached the side chains from the backbone skeleton, an approach previously 
adopted by Vishnyakov^i, and also cut the backbone skeleton into segments^^. In total the 
system then consists of A^^ sidechain segments and A*";, backbone segments. Taking into ac- 
count that the sulfonic acid groups are hydrophilic while the ethers and fluorocarbon groups 
are hydrophobic^^, we adopted the following general representation for a coarse-grained ar- 
chitecture of side chains and backbone segments. Each side chain contains = Np + 
monomers of which Np are hydrophilic and A'^^^ are hydrophobic. Its architecture is written 
as: Nnii+NpP. The backbone segments are fully hydrophobic with an architecture riH. In 
the majority of systems studied in this paper, the side chain architecture was 8H+2P and 
the backbone segment architecture was lOH. 

The simulations were then performed in several stages. The initial configuration was first 
prepared by using Monte-Carlo techniques to grow a mixture of polymer chain segments with 
fixed bending angle 6q and bond length Rq. After an initial equilibration the constraints on 6q 
and Ro were removed. Subsequently, MD runs were performed for time periods between 50ps 
and 500ps at constant volume V and constant temperature T until the system of polymer 
segments was fully equilibrated. The system temperature was controlled by coupling the 
system to a Langevin thermostat with a friction coefficient 7 = 0.1 and Gaussian white-noise 
force of strength QksT'j. The equations of motion were integrated using the velocity Verlet 
algorithm with a time step of 0.5 fs. We imposed standard periodic boundary conditions to 
our system, filling space with translational replications of a fundamental cell. Long-range 
electrostatic interactions were calculated using the standard Lekner summation algorithm^^. 

To check whether equilibrium had been reached we examined the proton pair distribution 
function g{r) in every picosecond time interval, and continued the run until no significant 
further change was observed. To avoid mistaking the trapping of the simulated system into 
a metastable glassy state for true equilibrium, we repeated each run with several different 
initial configurations. For each computer experiment we typically gathered system statistics 
over time intervals of a few nanoseconds. We believe this to be sufficient to allow for the 
self-organization of the side chains and backbones into clustered structures. 
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In the next stage of the simulations the segments were reassembled into the branched chain 
characterizing the original Nafion. This was achieved by the simultaneous introduction of a 
specific mutually attractive force between the ends of each backbone segment, which united 
them into a single chain, and a similar force acting between the tail monomer of the detached 
side chains and the fifth monomer of every backbone segment. This was applied in the form 
of a LJ potential of strength Slj — GhsT. To avoid the formation of star-like branched 
polymers, only a single occupancy of the backbone attachment sites was permitted. The 
simulation was then resumed, and run until a new equilibrium was attained. It was found 
that the bonds on the ends of all the backbone segments were occupied, and that a large 
majority (typically around 95%) of the side chains were reattached to a backbone segment. 
This fraction was a function of the ratio ip of the number Ng of side chains (and hence 
of sulfonate groups) to the number n x A^;, of backbone monomers. The quantity -i/; thus 
measures the concentration of sulfonate groups relative to the concentration of backbone 
monomers. The few side-chains that did not bond to the special sites on the backbone 
attached themselves more loosely to the backbone through the weaker regular hydrophobic- 
hydrophobic LJ interaction potential. 

We also introduce the absolute molar concentration of sulfonate groups, 77 = iVs/iVoV"ceii, 
with A^o being Avogadro's number and Vceii the volume of the simulation cell, as an 
alternative parameter with which to describe the sulfonate concentration, independently of 
the amount of backbone polymer. While it is true that in real dry Nafion, no significant 
free volume exists, and so the total number of monomers per unit volume is roughly fixed 
at the time of synthesis, in a simulation, however, it is possible to change the density of 
the material. It is of interest in understanding the mobility of ions in Nafion to know what 
structures would form in the presence of imposed free volume, were such a reduction in 
density of the material possible. 

In order to observe the effects of the dipolar nature of the sulfonate-proton group on 
the structural characteristics of side chain aggregates, we calculated the normalized pair 
distribution function g{r) of head group protons (which are the end group monomers of the 
side chains). The function g{r) indicates the relative probability of finding two protons at a 
separation distance r averaged over the equilibrium trajectory of the simulated system, so 
that 

.M = (Wiv.)^. (6) 
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Here dn+ is the number of protons located in a shell of thickness dr at the distance r from 
a fixed proton, and is the total number of protons in the system of volume V^g}}- 

In addition to this spherically averaged quantity, it is also useful to note some indications 
of the local anisotropy in order to distinguish between the Gierke and Gebel models. To 
this end we also evaluated the partial structure factors 
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Here j denotes the x, y and z components of the vector q, with qj — 2'!rnj/L, and the sums 
can proceed over the entire sample or over smaller volumes. 



IV. SIMULATION RESULTS 



The goal of these extensive coarse-grained molecular-dynamics simulations was to in- 
vestigate the relative importance of the various model parameters on the aggregation of 
the sulfonate (acidic) head groups. To achieve this we performed a series of simulation 
runs in which the key parameters were varied as summarized in Table H. These parameters 
included the absolute concentration and the relative concentrations of side chains and back- 
bone polymer, as well as the strength of the electrostatic interactions. This last quantity 
is characterized by a Coulomb coupling parameter defined as F = e^/ekjjTa = Xb/o'. Here 
Xb, which is known as the Bjerrum length, is about 56 nm for a dielectric constant e = 1 and 
T = 300 K. The coupling parameter F was varied between 2 and 160 by varying the system 
dielectric permeability e between 80 and 1, with higher F representing stronger electrostatic 
correlations in the system. In a typical simulation, a few thousand side chains are confined in 
a cubic box whose side is of length L — 30a, with a — 0.3 nm being the monomer diameter. 

Each simulation run was started from a randomly generated assembly of side chains and 
backbone segments within the system volume. However, in order to avoid packing problems, 
and to make this process effective for high monomer densities, a preliminary simulation of 
the same number of monomers was first performed in a larger box of side L' > L, and 
this was then followed by a gradual reduction of the system volume until the desired molar 
concentration was achieved. 
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Table I: List of key variables. 

Dj. dipole moment reduction factor 

^0 equilibrium bending angle for backbone and sidechain polymers 

Rq equilibrium bond length for backbone and sidechain polymers 

ke bending force constant 

kij stretching force constant 

kB Boltzmann constant 

T system temperature 

a monomer diameter 

£lj Lennard-Johns interaction parameter between monomers 

Qe normalized charge of sulfanates and ptorons 

e dielectric constant of medium 

Ng number of side chains 

Ni, number of backbone chain segments 

Njji number of monomers per side chain 

Nh number of hydrophobic monomers per side chain 

Np number of hydrophilic monomers per side chain 

rift number of hydrophobic monomers per backbone chain segments 
■0 ratio of the number of sulfonate groups to the number of backbone monomers. 

rj molar concentration of head groups 

A^o Avogadro's number 

L length of simulation box 

g{r) proton-proton pair distribution function 

77,+ local number of protons 

A^_l_ number of protons in the system 

Vceii volume of simulation box 

S{qj), j = x,y,z partial structure factors 

r coupling parameter in ionomer 

Xb Bjerrum length 
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Table II: Parameters used in combinations for simulation runs. Here r] is the molar sulfonate-ion 
concentration, Nh the number of backbone chain segments, L/a the simulation cell size in units of 
a = 0.35 nm, is the number of monomers per sidechain, ■0 is ratio of the number of sulfonate 



groups to the number 


of backbone 


monomers, 


and e is 


the dielectric 


constant of the medium. 


Runs 


ry(mol/l) 


Nb 


L/a 


Nm 




e 




Run 1 


0.8 





30 


10 





1, 15, 80 


1 


Run 2 


0.8 


500 


30 


10 


1 


1, 15, 80 


1 


Run 3 


0.8 


1000 


30 


10 


2 


1 


1 


Run 4 


0.8 


1500 


30 


10 


3 


1 


1 


Run 5 


0.8 


2000 


30 


10 


4 


1 


1 


Run 6 


0.8 


2500 


30 


10 


5 


1 


1 


Run 7 


1.6 





30 


10 





1, 15, 80 


0, 0.5, 1 


Run 8 


1.6 


iOOO 


on 
30 


iO 


i 


1 


0, 0.25, 0.5, 1 


Run 9 


1.6 


2000 


30 


10 


2 


1 


0, 0.5, 1 


Run 10 


0.2 


1000 


60 


10 


1 


1 


1 


Run 11 


0.5 


1000 


45 


10 


1 


1 


1 


Run 12 


2.5 


1000 


25 


10 


1 


1 


1 


Run 13 


0.16-4.8 





30 


10 





1 


1 


Run 14 


1.6 





30 


2-10 





1 


1 



A. The role of the strength of the electrostatic interactions 

In this first simulation the constraints imposed by the attachment of the side chains 
to the backbone were removed in order to concentrate on the effects of the electrostatic 
interactions. We then analyzed the proton-proton pair distribution function g{r) calculated 
for an assembly of unattached side chains for three different dielectric constants e of the 
ionomer medium, and for two different molar concentrations. The higher concentration, for 
which 7] = 1.6 mol/1, is characteristic of typical dry Nafion samples. Figure [T] shows the 
proton-proton pair distribution function for the parameters described in Table II as Run 1 
and Run 7. The first observation is that, as expected, small e, or, equivalently, high F, 
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facilitates clustering of the head groups, as shown in the red dot-dashed lines in Figure 
[U The graph of g{r) shows well-defined maxima at r/a = 1.5,2.0, and 3.4. A detailed 
investigation shows that this is a consequence of the formation of three simple multiplet 
structures: compact and linear multiplets, each created by two head groups, and a branched 
multiplet composed of three head groups, as shown in Fig. [21 

For e = 80 there is no correlation between head groups on different side chains except for 
steric hindrances, as seen in the continuous blue lines in Figured! A more surprising result is 
the effective repulsion between acidic head groups at intermediate dielectric permittivities, 
as seen in the green dashed lines in Fig. [H It appears that, in the absence of strong 
correlations, each head-group dipole prefers to stay in the bulk of the sample in order 
to increase its polarization energy. This is similar to the well known salt- depletion effect 
in colloidal physics, when salt pairs are effectively repelled from the surface of a neutral 
coUoid^^. 

The other unexpected result is the suppression of the correlations between the end groups 
at high molar sulfonate concentrations rj, as seen in the dot-dashed line and symbols in Figure 
[U^b). This can be traced back either to the reduction of the Debye screening length, which 
effectively shields the dipoles from each other, or to the enhancement of the hydrophobic- 
hydrophobic attraction in the system. The latter emerges from the hydrophobic tails of the 
side chain segments. 

To decide which of these two possible mechanisms was responsible, we added into the 
system described in Fig. lU^b) a small number of backbone segments in order to boost the 
effective hydrophobicity (or association energy) in our model system. One backbone segment 
of 10 monomers was added for each side chain unit, making ip = 0.1. After equilibration, the 
side chains were then loosely attached to the backbone segments and the backbone segments 
were mostly united in the manner previously described, whereupon the system was further 
equilibrated. The results, shown as Run 2 in Figure [H clearly indicate that the backbone 
segments have a minor impact on the head group distribution function. The system has 
almost the same microphase structure as observed for the original side chain system. A 
fivefold increase in the backbone concentration in Figure [3]^b) has only a small effect on 
g{r). Therefore, it is the dipole-dipole interaction that tends to destroy the correlations 
when the molar fraction of sulfonates is high. This effect takes place when the screening 
length becomes less than the mean separation distance between the dipoles. This has some 
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relevance to the search for ionomer membranes with superior transport properties: whereas 
for better proton conductance one needs a higher sulfonate molar fraction, the latter hinders 
sulfonate aggregation, which is vital for water-sulfur channel formation in wet ionomers. 

Typical snapshots of the simulation box, shown in Figure HI make visible the process of 
microphase separation in ionomer systems. The tips of the side chains, which correspond 
to end group protons, are shown as balls, whereas the rest of the side chains are plotted as 
lines. The backbone segments, drawn also as lines, mostly cluster together to create large 
aggregates. The side chains either anchor onto the surface of the backbone aggregate, or 
create their own small clusters. 

Figure presents a 3-dimensional density distribution of monomers for Run 5. The areas 
of low (high) density colored in blue (red) show the channel structure of the aggregates. The 
polymer matrix is composed of hydrophobic tails of side chains and backbone segments. The 
hydrophilic head groups form the walls of these channels. Though this 3D density picture 
differs from the long-bundle structural model of Nation^, it supports the model of a dynamic 
system of non-spherical clusters with temporary inter-cluster bridges^. 

B. Structure dependence on sulfonate molar concentration rj. 

In the previous section, in which we varied the strength of the Coulomb interactions, we 
observed that in the systems with the higher sulfonate concentration, the Coulomb correla- 
tions were suppressed. This suggested that a study of the morphology of aggregates over a 
wider range of concentrations would be helpful to understand the physics behind sulfonate 
cluster formation. Simulation results for g{r) for systems with sulfonate concentrations r) 
varying between 1.6 and 2.5 mol/1 are shown in Fig. [6]for Runs 8 and 10-12. Results for the 
much wider range of 0.16 to 4.8 mol/1 are shown in Fig. U\ but for these it was necessary to 
delete the backbone segments entirely in order to accommodate the very large concentration 
of side chains at the upper part of the range constituting Run 13. 

The dipole-dipole correlations are seen to depend strongly on the molar concentration of 
acidic groups. For a typical dilute case. Run 10 in Fig.[6l the first peak of g{r) reaches a value 
g{r = 1.5a) = 37, indicating a strong condensation of nearly all the side chains. A possible 
explanation for this effect can be drawn from the clustering features of pure dipolar liquids. 
In a dilute system, where the average interchain distance is larger than the gyration radius Rg 
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of the side chain, the LJ interaction between the neighboring chains is neghgible. The main 
contribution to the energy in this case comes from the electrostatic dipole-dipole interaction 
of head groups. Although this interaction decays with a distance dependence, it is not 
weak over distances comparable with the simulation box length L, and is capable of forcing 
the chains to aggregate. However, for a dense system, where the LJ contribution to the 
energy becomes comparable to the dipolar term, the side chains start to lose their rotational 
entropy, and thus cannot reorient themselves in order to minimize the electrostatic energy 
of the aggregate. While the heights of the maxima of g{r) gradually decrease as the system 
becomes more tightly packed, their positions do not change and remain fixed at r = 1.5a, 2a, 
and 3.4(7. 

The line with circles in Fig. [7] corresponds to the case where the concentration of back- 
bone segments was double the case for Run 8, and corresponds to Run 9. It is seen that the 
proton pair distribution function is hardly affected by the existence of the extra hydropho- 
bic segments. However this conclusion is not necessarily true for a system with the more 
attractive LJ parameters for HH and HP interactions that will be considered in section D. 
Both the snapshot picture in Fig. [8] and the 3D monomer density in Fig. [9] show microphase 
separation. 

Another way to examine the head group clustering is to calculate the partial structure 
factors S{qj) {j = x, y, z) of the protons. Because these did not exhibit any marked 
anisotropy, an averaged S{q) is plotted in Fig. [101 The position of the peak at qa = 4.2 
is associated with a correlation distance, d = 1.5a in the sample structure via the Bragg 
relationship d = 2TT/q. A sharp growth at small q of the values of S{q) in Fig. [TO] results 
from microscopic aggregation of hydrophilic and hydrophobic monomers in the system. 

To summarize the tendencies of dipole-dipole correlations presented in this and the pre- 
vious section, we conclude that the hydrophilic head group aggregation takes place prefer- 
entially at low dielectric permittivity e and low densities of sulfonate groups. 

C. Role of side chain size Nm- 

Although the bulk of this work describes the conformational structures of side chains of 
length Nm=10, it is worth while investigating systems with shorter side chains, since recent 
experimental show a high proton transfer rate in ionomer membranes with shorter pendant 
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segments. 

This section contains simulation results for different side chain lengths, A^^ = Np + Nh- 
The length of the hydrophilic part (the number of hydrophilic monomers) is fixed at Np = 2 
and the length of the hydrophobic tail (the number of hydrophobic monomers) is varied 
between and 8 (which corresponds to the variation of total side chain length between 2 
and 10) in order to reveal how the side chains' hydrophobic tail affects the pair distribution 
function of protons. No backbone material was included in these systems. Simulation 
results for Run 14 are shown in Fig. [TTl Note that the case = 2 corresponds to a system 
of dipolar dumbbells. Even a short hydrophobic tail with Nm = 4 strongly suppresses 
the aggregation of head groups. Thus the hydrophobicity of side chains, or the length 
of hydrophobic tail, is a strong obstacle to head group clustering. In other words, the 
dipole-dipole attraction between sulfonates of longer side chains is not sufficient to produce 
a nanophase separation. This conclusion is in contrast to the suggestion of Eisenberg^ 
that each ionic group in the vicinity of a cluster of sulfonate groups is drawn in by dipolar 
attraction to form a depleted zone free of charge. 

D. lonomer morphology for strong hydrophilic and hydrophobic interactions. 

One intriguing question in fuel cell technology is how the morphology of nanophase sep- 
aration in PEM membranes depends on the polymeric structure of Nation. In other words, 
what is the optimal comb-like structure for the ionomer to be a good material for a fuel cell? 
Using our segmental model for Nation we can mimic different comb-like Nation molecular 
conformations by changing the parameter ip, which denotes the ratio of the number of sul- 
fonate groups to the number of backbone monomers. The impact of on the distribution 
function g{r) was shown to be negligible for Elj = SksT for HH and Elj = ^kpT for both 
HP and PP interactions in Fig. [3] for Runs 1-6. In order to explore more generally the effect 
of changing ip, we exaggerate the hydrophilic and hydrophobic interactions by tripling the 
interaction constants Elj. The Results for Runs 1-3 presented in Fig. [12] for e^j = OksT for 
HH, and e^j = SksT both for HP and PP interactions, now show that the addition of back- 
bone chain segments brings more conformational disorder into the head group correlations. 
This is also evident from the system snapshots for Run 3 in Fig. [12] for the two different 
strengths of Lennard- Jones interaction parameters. In other words, a strong interaction 
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between hydrophobic and hydrophihc monomers can disrupt the side chain aggregation, 
causing all the side chains to be evenly distributed on the surface of backbone aggregations, 
as in Fig. [T51(b). This distribution will have a globular, lamellar, or cylindrical-rod form. 

E. End group attraction between side chains 

To complete this study of the effects of interactions on morphology, we look at the effect of 
both decreasing and increasing the interaction between head groups. We have already seen in 
section IIIA some of the consequences of reducing the Coulomb interaction by increasing the 
dielectric constant e. In the absence of interactions other than the dipole-dipole interaction, 
this is equivalent to a reduction in the strength of the dipole moments. In Fig. [H] we show 
the effect on the structure factor S{q) of the protons of a reduction in the strength of the 
dipole moments by various factors Dr (Run 8). When the parameter Dr is small, there is 
no phase separation in the ionomers, and the function S{q) has a liquid-like structure with 
a single maximum. The behavior of S{q) is changed for an unreduced dipole strength Dr=l 
in a manner similar to the behaviour of S{q) in the less dense system shown in Fig. [TOl 
there is a strong peak at qa=A.2, and a sharp growth at small qa. The latter, again, is an 
indication of a strong microscopic separation of hydrophilic and hydrophobic monomers in 
the system. 

It is clear that the effective attraction between the acidic groups of side chains can- 
not be accounted for by the sulfonate-hydrogen dipole moment alone. The sulfonate end 
group has its own internal dipole moment due to the charge distribution between the sul- 
fur atom and the three oxygen atoms bound to it. Partial charges on the other side chain 
monomers can also contribute to the effective attraction between head groups. To take all 
these contributions into account at a phenomenological level we introduce an artificial at- 
tractive interaction between the end groups of the side chains, in an approach similar to that 
employed by Shirvanyana^. For this purpose we use a Yukawa-type attractive potential, 

F(r) = — (Z^e^/er) exp (— r/r/)) for r > a ; F(r) = oo for r < a. 

Here Z is an effective Yukawa charge for the head group sulfonates. For the sake of simplicity, 
we restrict ourselves to the case of constant Debye screening length r^y = Aa. In fact, the 
parameters Z and determine the stability, size and shape of the aggregates forming in the 
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system. The simple attractive Yukawa term can serve as a convenient basis for interpreting 
all the important features that can be observed in real systems or predicted by statistical 
mechanical theories. 

In the absence of electrostatic interactions in the system, i.e. when the head group dipole 
moment vanishes and = 0, the ionomer system becomes spatially inhomogeneous. The 
pair distribution functions g{r) for Runs 7-9, plotted in Fig. [15] for Z = 0.1 (case (a)) and 
Z = 0.5 (case (b)) show strong structuring. The aggregation is substantially stronger for 
a system missing backbone polymer (Run 7) than for a system with a full complement of 
backbone monomers (Run 9), as seen from a comparison of the respective lines in Fig. [151 
There is a clear tendency to micellar arrangement, which is characterized by some short- 
range order. Snapshot pictures, given in Fig. [T^l show formation of a space-filling web-like 
network of cluster-like aggregates surrounded by neutral polymer sections. The size of these 
clusters grows as the density of side chains increases. A comparison of Figs. [TBT a) and [TBT b) 
reveals that the effective Yukawa charge Z, or in other words, the strength of the the screened 
Coulomb attraction, has a strong impact on the shape of the clusters. The structure of the 
system for a moderate Yukawa attraction of Z = 0.1 between the end groups is mobile: the 
aggregates can diffuse, decay and reappear again. They can also create temporary bridge-like 
connections. Therefore the nanophase separation of the head groups, shown in Fig. [TBT a) 
can be viewed as an example for a Gierke-like sulfonate aggregation. However, a strong 
attraction between the end groups of the side chains, the case Z = 0.5 in Fig. [T6](b). results 
in a frozen network of globular clusters of diameter roughly 1.5 nm in the polymer matrix. 
This cluster picture resembles the ion segregation phenomenon proposed by Eisenberg^^. 

When the dipole moment of the sulfonate head groups is restored, in this case with a 
value of Dr = 0.5, the globular micelles form super-clusters of diameter 3-4 nm , as seen 
is Fig. [TTT b) for Z = 0.5. On the other hand, this nonzero dipole moment facilitates the 
formation of percolated cluster structures for Z = 0.1, as in Fig. [TTT a). The inclusion of 
the dipole moment of the sulfonate head groups thus changes the conformational structure 
of clusters from a Gierke model of spherical aggregations to a Gebel model of elongated 
aggregations. A snapshot picture of the same system but with no artificial attraction be- 
tween the head groups is given in the Fig. [181 It is clear that the introduction of a screened 
Coulomb attraction results in the formation of compact and interconnected clusters of sul- 
fonate head groups. 
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V. CONCLUSIONS 



In this study we have explored the issue of nanophase separation in dry Nafion® films 
by means of simulations within a simple dipolar-head-group model for the side chains. Our 
procedure involved detaching the side chains from the backbone and cutting the backbone 
into segments to speed equilibration, and then reassembling the macromolecule by means of 
a strong imposed attractive force between the cut ends of the backbone, and between the 
non-ionic ends of the side chains and the midpoints of the backbone segments. 

As we vary the interaction parameters and dielectric permittivity of the medium we 
find different structures to be formed by the head groups. Small permittivities facilitate the 
clustering of head groups with a well defined maximum in the proton-proton pair distribution 
function g{r) that indicates three simple multiplet structures: compact and linear multiplets 
each created by two head groups, and a branched cluster composed of three head groups. At 
higher permittivity no correlations between head groups are obtained. An effective repulsion 
between acidic head groups at intermediate dielectric permittivities was revealed from our 
data. We assume that, in the absence of strong correlations, each head-group dipole prefers 
to stay in the bulk of the sample in order to increase its polarization energy. 

For the ionomer interaction parameters considered in this work the microphase structure 
of the head groups was found to have little dependence on the concentration of backbone 
segments. Since the molar fraction of sulfonates in a fuel cell membrane is known to be a 
crucial parameter for its performance, this factor was also studied. We find that for a dilute 
case a strong condensation of nearly all side chains takes place. The main contribution to the 
Hamiltonian in this case comes from a long-range dipole-dipole interaction of head groups 
which forces the chains to aggregate. However, for a dense system, where the hydrophobic LJ 
contribution to the Hamiltonian becomes comparable with the dipolar term, the side chains 
cannot reorient themselves in order to minimize the electrostatic energy of the aggregate. 
While on the one hand, for high proton conductance one needs higher sulfonate molar 
fractions rj, on the other hand we find that a high value of this parameter rj hinders the 
sulfonate aggregation that is vital for water-sulfonate channel formation in ionomers. 

Overall, our simulations show that the dipole-dipole interaction is the main driving force 
for the aggregation of sulfonates. When the dipole moment of sulfonates is small, so that 
Dr ^1, there is no phase separation in the ionomers. In contrast, for = 1 a microscopic 
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aggregation of hydrophilic and hydrophobic monomers takes place in the system. 

The most significant feature learned from the simulations is the fact that the dipolar- 
head-group approximation is not sufficient to drive the compact side chain cluster forma- 
tion suggested by experimental measurements. An additional artificial attraction between 
the sulfonates, introduced as a screened Yukawa interaction potential, boosts formation of 
nanometer-sized clusters dispersed throughout the system volume. The final size of these 
clusters, and the geometry of the cluster network they create, depend on the dipole mo- 
ment of the sulfonate head groups for fixed Yukawa parameters. As the strength of 
increases from zero to one, the conformational structure of clusters changes from a Gierke 
model of spherical aggregations to a Gebel model of elongated aggregations. 

More work is required to develop a deeper understanding of the mechanisms behind the 
aggregation picture observed experimentally. Because the shape and size distribution of 
the aqueous pores in wet membranes is not completely known, it will be useful to extend 
our preliminary results for the dry Nafion considered in this paper to systems with water 
present and with a more detailed chain structure. In future work, we plan to extend the 
model presented here to more complex models. In particular we will employ partial charges 
and Lennard- Jones parametersi^ii^>2^ for side chain and backbone monomers of Nafion-like 
oligomers. The issues of how the partial charges of ionomer molecules, in fixed and free 
proton models for side chains, affect the nanophase structure of the PEM materials will be 
addressed. A complete understanding of the role of head group charges on the cluster confor- 
mations of side chains can provide helpful guidelines for understanding existing membranes 
and designing new membrane materials for fuel cells with improved characteristics. 

Acknowledgments 

This work was supported by the US Department of Energy under grant DE-FG02- 
05ER46244, and was made possible by use of facilities at the Case ITS High Performance 
Computing Cluster and the Ohio Supercomputing Center. 



^ M. I. Perry, T. F. Fuller, J. Electrochemical Society 149, S59 (2002). 

2 T. D. Gierke, G. E. Munn, F. C. Wilson, J. Polym. Sci. Phys. Ed. 19, 1687-1704 (1981). 



19 



^ W. Y. Hsu, T. D. Gierke, Macromolecules 15, 101-105 (1982); ibid J. Memb. Sci. 13, 307-326 
(1982). 

^ S. J. Paddison, Annu. Rev. Mater. Res. 33, 289-319 (2003). 

5 G. Gebel, O. Diat, Fuel Cells 5, 261 (2005); 1. Rubatat, O. Diat, G. Gebel, Physical chemistry 
and Soft Matter 106, 1 (2004). 

6 S. J. Paddison, J. A. Elliott, J. Phys. Chem. 109, 7583 (2005). 

^ M. Eikerling, S. J. Paddison, L. R. Pratt, T.A. Zawodzinski Jr., Chem. Phys. Letters 368, 108 
(2003). 

8 K. D. Krcucr, S. J. Paddison, E. Spohr, M. Schuster, Chem. Rev. 104, 4637 (2004). 

9 T. Li, A. Wlaschin, P. B. Balbuena, Ind. Eng. Chem. Res. 40, 4789-4800 (2001). 
^° K. A. Mauritz, C. E. Rogers, Macromolecules 18, 483 (1985) 

M. Eikerling, A. A. Kornyshev, U. Stimming, J. Phys. Chem. B 101, 10807 (1997). 
^2 A. S. loselevich, A. A. Kornyshev, J. H. G. Steinke, J. Phys. Chem. B 108, 11953 (2004). 
^3 M. Laporta, M. Pegoraro, L. Zanderighi, Phys. Chem. Chem. Phys. 1, 4619-4628 (1999). 

P. G. Khalatur, S. K. Talitskikh, A. R. Khokhlov, Macromol. Theor. Simul. 11, 566-586 (2002). 

15 T. A. Zawodzinski Jr., M. Neeman, L. O. Sillerud, S. Gottesfeld, J. Phys. Chem. 95, 6040-6044 
(1991). 

16 A. Vishnyakov, A. V. Ncimark, J. Phys. Chem. B 105, 7830 (2001); ibid 9586-9594 (2001). 
1^ J. T. Wescott, Y. Qi, L. Subramanian, T. W. Capehart, J. Chem. Phys. 124, 134702 (2006). 
1^ S. Yamamoto, R. Jinnouchi, Sh. Yamakawa, Sh. Hyodo, 14th International Conference on the 

Properties of Water and Steam, Kyoto, Japan, p.411 (2004). 
1^ E. Spohr, Molecular simulations 30, 107 (2004). 

20 S. J. Paddison, T. A. Zawodzinski, Solid State Ionics 113-115, 333-340 (1998). 

21 A. Vishnyakov, A. V. Neimark, J Phys. Chem. B 104, 4471 (2000); D. Rivin, G. Meermeier, N. 
S. Schneider, A. Vishnyakov, and A. V. Neimark, J. Phys. Chem. B 108, 8900 (2004). 

22 S .C. Glotzer, W. Paul, Annu. Rev. Mater. Res. 32, 401-436 (2002). 

23 M. Mazars, J. Chem. Phys. 115, 2955-2965 (2001). 

2^ E. AUahyarov, J. -P. Hansen, A. Lous, H. Lowcn, Europhys. Letters (2005). 

25 V. Arcella, A. Ghielmi, G. Tommasi, Ann. N. Y. Acad. Sci. 984, 226 (2003). 

26 A. Eisenberg, Macromolecules 30, 7914 (1997). 

2^ D. G. Shirvanyanz, A. S. Pavlov, P. G. Khalatur, J. Chem. Phys. 112, 11069 (2000). 

20 



S. S. Jang, V. Molinero, T. Cagin, W. A. Goddard III, J. Phys. Chem. B 108, 3149-3157 (2004). 



21 



Figure 1: (Color online) Proton-proton pair distribution function g{r) for (a) Run 1 (low density), 
and (b) Run 7 (high density) for e=l (dot-dashed line), 15 (dashed line) and 80 (full line). Simu- 
lation results for e=l are scaled down by a factor of four to fit the figure dimensions. The e = 1(*) 
(symbols) in (b), given for comparison, corresponds to the e = 1 line in (a). 
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Figure 2: (Color online) Simple multiplet structures in side chain aggregates. Upper row: compact 
(left ) and linear (right) multiplets. Bottom row: branched multiplet. The hydrophobic tails are 
shown only for a compact multiplet. The dashed objects C and D for the branched cluster indicate 
another possible configuration of objects A and B respectively. Positively charged hydrogen atoms 
are shown as yellow circles. 
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Figure 3: (Color online) (a) Proton-proton pair distribution function g{r) for Run 2 and different 
dielectric permittivities. e=l (dot-dashed line), 15 (dashed line) and 80 (full line). Results for 
Run 1, not shown here, are very close to the results for e=l. Simulation results for e=l are scaled 
down by a factor of four to fit the figure dimensions, (b) Proton-proton pair distribution function 
g{r) for different amounts of backbone polymer (Runs 2 and 6). Results for Run 3-5 (not shown 
here) that correspond to intermediate amount of backbone polymer, lay between these two curves. 




Figure 4: (Color online) Typical snapshots of a system for Run 4 (case (a)), and Run 6 (case (b)). 
Spheres represent the end-groups (protons) of side chains. The sulfonate group and neutral chain 
sections are drawn by lines (in red in online version). Different bead colors correspond to different 
bead altitudes, with a blue color for low-altitude beads (at the bottom of simulation box) and a red 
color for high-altitude beads (at the top of simulation box). The size of all the structural elements 
is schematic rather than space filling. 





Figure 5: (Color online) A short-time average of the 3D density of monomers for Run 5. The hy- 
drophobic tails of side chains and backbone segments create a solid matrix, whereas the hydrophilic 
head groups forms the walls of the channels. 




Figure 6: (Color online) Proton-proton pair distribution function g(r) for different box sizes L, 
and hence sulfonate concentrations ij, for Runs 8 and 10-12. The first peak of g{r) for Run 10 (full 
line, blue in colored figure) has a value g{r = 1.5a) = 37, out of the scale of the plot. 
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Figure 7: (Color online) Proton-proton pair distribution function g{r) for different sulfonate molar 
concentrations rj, Run 13. (a): 0.16 mol/1 < rj < 1.6 mol/1. (b): 1.6 mol/1 < i] < 4.8 mol/1. The 
case rj = 1.6(*) (line with circles in (b)) is for Run 9. The first peak of g{r) for rj = 0.16 (rj = 0.32) 
mol/1 (full blue line (dashed green line) in (a)) has a value g{r = 1.5a) = 53 {g{r = 1.5a) = 22), 
which is out of the scale of the plot. 
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Figure 8: (Color online) Typical snapshots of system for Run 13 containing rj = 0.4mol/l (case (a)) 
and 7] = 3.2 mol/1 (case (b)) head groups. 




Figure 9: (Color online) A short-time average of 3D density of ionomer monomers for Run 13 and 
rj = 2.4 mol/1. The areas of low density are cut away. 
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Figure 10: (Color online) Structure factor S{q) of the head groups for Run 3 versus wave number 
qa. 
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Figure 11: (Color online) Proton-proton pair distribution function g{r) for a system of side chains 
for Run 14. The number of monomers per side chain segment is changed from 2 to 10. 




Figure 12: (Color online) Proton-proton pair distribution function g(r) for Run 1 (solid line), Run 2 
(dashed line) and Run 3 (dot-dashed line). Thick lines: elj = Qfc^T for hydrophobic-hydrophobic 
(HH) interaction, elj = SksT for hydrophobic-hydrophilic (HP) and hydrophilic-hydrophilic (PP) 
interactions. Thin lines: elj = SksT for HH, elj = IksT for HP and PP interactions. 
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Figure 13: (Color online) Snapshot pictures of systems for Run 3 and n^V =2- Case (a): slj = 
SksT for HH, elj = ^ksT for HP and PP interactions. Case (b): e^j = ^ksT for HH and 
elj = SksT for HP and PP interactions. In the case (b) the sulfonate groups are mostly on the 
surface of clusters. 





Figure 14: (Color online) Proton-proton structure factor S{q) for Run 8. Dipole strength Dr is 
scaled down from 1 to 0. 
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Figure 15: (Color online) Proton-proton pair distribution function g{r) for Runs 7-9 and Z = 0.1 
(a) and Z = 0.5 (b). Z is a fictitious charge on the head groups which causes an artificial Yukawa 
like attraction between sulfonates. Thick (thin) lines are for Dr = {Dr = 0.5). Note that the 
scales on y-axis for short-range and long-range behavior of g{r) (shown as an inset) are different. 
For other details see the text. 
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Figure 16: (Color online) Snapshot pictures of systems for Run 9, rihip =2, Dr=0, Z=0.1 (case (a)) 
and Z=0.5 (case (b)). 





Figure 17: (Color online) Snapshot pictures of systems for Run 9, rihtp =2, Dr=0.5, Z=0.1 (case 
(a)) and Z=0.5 (case (b)). 





Figure 18: (Color online) Snapshot of a system for Run 9, rihip =2, D=0.5, Z=0. 
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